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(54) Modelling an underground reservoir 

(57) A model simulating the behaviour of an underground reservoir (eg of oil) is defined on the basis of 
initially known geological information and available data. It enables production forecasts to be made using 
parameters representative of the initial data. The method essentially consists in defining one or several 
possible scenarios of developments during production by creating, for each of these scenarios, fictitious new 
production data corresponding to hypotheses about the future states of the reservoir. For each scenario 
considered, it is checked whether the available parameters of the simulation model can be adjusted, given the 
constraints of the initial geological model so that the simulation model reproduces both the measured 
production data and the added data. The method may be used to quantify the uncertainties regarding 
production forecasts by finding the min/max extremes of future production values. 
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METHOD FOR PREDICTING BY AN INVERSION TECHNIQUE 
CHANGES THAT WILL OCCUR DURING PRODUCTION IN 
AN UNDERGROUND RESERVOIR 

The present invention relates to a method for 
predicting changes that will occur during production in 
an underground reservoir, and in a reservoir containing 
hydrocarbons in particular, by means of a scenario 
inversion technique. 

BACKGROUND OF THE INVENTION 

Inversion techniques are used extensively in the 
field of reservoir engineering. The following fields of 
application may be cited as examples: 

- At the laboratory stage, they are used to determine 
various parameters representing the behaviour of rocks 
with respect to fluids. These parameters may be, for 
example: the absolute and relative permeability, 
capillarity curves, etc. . 

- Inversion techniques are also currently used to 
interpret well tests. In this case, the inversion 
parameters are, for example, the permeability of one or 
several facies, the geometric limits of a geological 
structure, coefficients pertaining to the productivity 
index of a well, etc.. 

- In studying a reservoir, inversion is used to match 
the response of a numerical simulator with available 
production measurements (or the "production history")- 
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The parameters may be, for example, the porosity of the 
rocks, their absolute and relative permeabilities, the 
productivity indices of the wells, etc.: 

5 Technique of matching the production history of a 

reservoir 

Inversion techniques represent a very active field of 
research in reservoir simulation- The process of matching 

10 a production history by adjusting parameters of a 

simulation model is also referred to by specialists as 
"history matching" . It mainly consists in finding a set 
of parameters governing the fluid flow equations which, 
when integrated in a numerical simulator, enable observed 

15 or indicated data to be found again. 

By way of example, the data observed in each well are 
mainly pressure, the composition of the fluids or the 
flow rate of the different phases. Any combination of 
these variables may also be used and in particular the 

20 water cut, the G.O.R. (Gas Oil Ratio), etc. Generally 

speaking, there may be many solutions to the inverse 
problem which will allow the production history of a 
reservoir to be reconstructed by a process of adjustment. 
It is therefore vital to incorporate what geological data 

25 is available in order to limit the range of possible 

solutions. By taking these initial geological data into 
account, the production forecasts produced by a 
simulation will be more realistic and more reliable. 

In order to construct the initial model, any 



information that is available is incorporated: raw or 
interpreted data, geological studies, seismic findings, 
etc- The physical data of a reservoir can be integrated 
by considering, for example: 

- the structure of the sedimentological units; 

- the limits of variations in petrophysical values 
(porosity, permeabilities, etc.) associated with litho- 
facies; and 

- statistical data pertaining to mean values, 
standard deviations, spatial correlation, etc. . 

The process known as history matching conventionally 
consists of the following steps: 

a) A simulation model is constructed using the 
initial geological data as a basis, integrating the 
maximum amount of data available . 

b) Significant parameters of the model are selected 
for the inversion process by considering: data relating 
to the behaviour of the well as a function of these 
parameters, their qualitative influence on production, 
the initial uncertainties associated with these 
parameters, etc.. 

c) Adjustments are made using the model parameters in 
order to reconstruct the production history or the 
observed data. Having established a set of parameter 
values, a direct simulation enables a comparison to be 
made between the predicted results and the observations. 
The most commonly used method is a process of trial and 
error in which the experience of the reservoir engineer 
comes into play: the parameter values are adjusted on the 



basis of the reservoir data and the understanding of its 
dynamic behaviour. 

It is sometimes possible to speed up this adjustment 
stage by using an automatic process that enables an 
iterative modification of the value of the parameters so 
as to obtain a better match between the result of the 
calculation and the observation data. 

The disadvantage of this conventional, largely 
empirical method is that the inversion process does not 
directly quantify the effect that the. initial knowledge 
of the reservoir has. In practice, it is possible to 
match a production history whilst deviating from the 
admissible physical quantities. Furthermore, 

uncertainties in the production forecasts are not 
directly quantified during the inversion process. 
Sensitivity studies are sometimes carried out on the 
parameter values obtained by the adjustment process in 
order to ascertain their effect on the forecasts but 
these studies are not an integral part of the inversion 
process . 

The process of history matching can be formalised as 
follows : 

Physical situations are modelled using a mathematical 
model designated as g, which in this case is the 
simulator of fluid flows in the reservoir. This model is 
dependent on ix, parameters (vector x) and time t. The 
physical situations modelled can then be written in the 
form d - g(Xft) . If d^ represents a set of n^ 
observations corresponding to the measurements of d at 



different times t if then the following conditions can be 
formulated so that the model reproduces the observations : 
d« » {g(x,t t ) ,i = 1,1^} 

In order to define a criterion for adapting the model 
to the observations, a function can be defined in the 
conventional manner of the "least squares", expressed, 
for example, as follows: 

F m (x) = l/2Si=l / n m { (di-gfx,^)) 2 } 

The parameters retained are those which minimise the 
function F m (x) by taking account of the constraints 
imposed on the parameters . 

It should be noted that the quality of the 
predictions produced by this procedure may sometimes be 
open to question for two main reasons. The first is that 
the initial or a priori knowledge of the reservoir is not 
quantified, which means that it may disappear during the 
minimisation processes. The second is that there may be 
several solutions which minimise the function. This 
being the case, in order to improve the quality of the 
predictions, it is important to take account of any 
information that is available about the geological model 
as well as any statistics that will characterise these 
parameters during the adjustment process. 

Bayesian technique of inversion of a simulation model 

In order to take account of this initial information, 
an inversion formalism known to specialists as "Bayesian 
Inversion" can be used. This is described, for example, 



by: 

Floris F.J.T. et al, in "Flow Constrained Reservoir 
Characterization using Bayesian Inversion" , 4th ECMOR, 
Norway, 7-10 June 1994. 

This formalism, which is known to specialists, allocs 
account to be taken of a priori data pertaining to tlie 
parameters in the form of probabilistic models, errors 
made in the measurements along with their relative 
probabilities as well as errors caused modelling flows 
using a numerical simulator. Given the calculating means 
currently in use, Bayesian formalism is unfortunately 
very costly in terms of calculating time and is not vexry 
practical where an entire field is under consideration 
since it requires the flow simulator to "sweep" all the 
possible parameters in order to evaluate the probability 
of each solution. 

The underlying principle of this Bayesian inversion 
as applied to flow simulation in a reservoir is to modify 
the initial geological model by using the production data 
in order to improve the quality of this initial model. 
The main steps are : 

a) Constructing an initial geological model. 

This is an "a priori" model of the probabilistic 
type. The reservoir is described on the basis of a set of 
parameters together with their associated probability 
density functions p(x) . Mean values for these a priori 
parameters {designated as a ^d standard deviations 

(designated as ax) can be used to model this probability 
density function, designated as pdf . With Gaussian type 



uncertainties, modelled by a covariance operator 
designated as Cx, the a priori function pdf can be 
expressed as : 

p(x) = Cte exp{-(Xa - x) T C/ 1 (x^ - x) } 

In order to express this a priori data in an 
objective function, a new term, designated as F x (x) , has 
to be added to the previously defined term F m (x), and can 
be expressed as follows in the case of Gaussian 
uncertainties : 

F x (x) :l/2{- (2^ - x) T C/ 1 (*a - x) 

b) Using the observed data: 

the observation data are used to calculate the a 
posteriori probability density function p(x/d=dm) in 
accordance with the Bayes rule: 

P(x/d=dj = (pCfl-d./*) .p(x))/p (d=cL) , 

where pf^d^/x) is the likelihood function which 
quantifies the differences between the observed values 
and the calculated values d for a given* set of parameters 
x. If the uncertainties about the measurements and about 
the model are Gaussian arid are modelled by a covariance 
operator C m , this likelihood function can be expressed as 
follows : 

Pte-i/*) - Cte exp{-(d.-d) T C B - l (d m -fl) } 
In order to obtain a statistical interpretation in 
terms of probability densities pdf, it is necessary in 
practice to determine the a posteriori function for the 
full range of parameters. If the maximum a posteriori 
probability only is considered, the parameters 
corresponding to this optimum are those which minimise 
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the objective function F=F m (x) + F x (x) . 

To sum up, for each set of parameters, the associated 
probability p{x) is increased if the response of the 
model g(X/t) provides values close to those observed c^, 
i.e. when the likelihood function pfd^d^/x) is high. 

The conventional method of history matching and the 
Bayesian inversion technique, mentioned and discussed 
above, are used to adapt a simulation model to actual 
measurable data. The major disadvantages of the first 
method are the possibility of losing a part of the 
initial information and a high degree of uncertainty as 
to the reliability of the predictions. Bayesian formalism 
allows the uncertainties pertaining to the geological 
parameters to be characterised more accurately and takes 
account of the a priori geological data. However, 
inversion of the Bayesian type is, in practice, very 
difficult to apply to real case studies, the cost 
incurred by the number of simulations needed being 
prohibitive . 

Furthermore, these methods are not intended as a 
response to the problem of quantifying uncertainties 
linked to production forecasts. Conventionally, the 
uncertainties are estimated a posteriori, after the 
inversion phase. Forecasting uncertainties, however, is 
crucial in drawing up a reservoir development programme, 
for which economic criteria based on the reliability of 
the forecast are required. 

The method of the invention allows the quantification 
of production uncertainties to be directly integrated in 
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the inversion process in a same procedure. The objective 
of this method is both to: 

- enable optimum matching of a simulation model to a 
production history by taking account of all the a priori 

5 information , and 

- enable uncertainties about production forecasts to 
be directly quantified and test hypotheses on the 
evolution of these production forecasts. 

The method enables a reservoir engineer to test 

10 production forecast scenarios and thus better validate 

his model against the initially known data of ttie 
reservoir. Consequently, it is possible to verify whether 
a production scenario is compatible with the geological 
data or quantify the production uncertainties inherent in 

15 the geological uncertainties. 

Forecasts of the likely evolution of production of an 
underground zone containing fluids, such as a hydrocarbon 
reservoir, can be produced by the method of the invention 
by carrying out the following steps : 

20 a) In a conventional initial step, already described 

in connection with the commonly used inversion methods, 
an initial geological model defined by a set of 
parameters is constructed by integrating all the 
available data. At this stage, conventional inversion 

25 methods can be ' used to determine the values of parameters 

which will allow the response of the model to be adjusted 
to the production history; 

b) One or several possible development scenarios for 
production are defined. For each of these scenarios, 



10 



fictitious data are added to the observed data (or 
"actual" data) at given instants in the future 
corresponding to these hypotheses ; 

c) For each scenario, the parameters of the said 
5 model are adjusted to reconstruct both the actual 

measurements derived from the production history and the 
data added at future times and an inversion procedure is 
applied in order to obtain a new set of parameters x for 
each scenario, which characterise a modified geological 
10 model and correspond to the hypotheses of the scenario; 

d) Direct simulations are then produced for each of 
these modified geological models in order to establish 
new production forecasts. 

It is possible to use conventional procedures in each. 

15 inversion phase. However, it is of advantage to use an 

automatic inversion procedure in which an objective 
function is minimised by incorporating the initial data 
associated with the geological model and to use the 
Bayesian formalism to take account of these initial data- 

20 To this end, the method may include a combination of 

an inversion formalism of the conventional Bayesian type 
and an optimisation algorithm based on the use of the 
gradients method, enabling the derivatives of the 
production forecasts to be obtained in relation to the 

2 5 parameters. By combining this method with a Bayesian 

formalism, it is possible to reduce the number of 
operations needed to obtain the optimum probability by a 
considerable degree, since the advantages of the Bayesian 
inversion (which takes account of the initial data) are 
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combined with the efficiency of inversion algorithms 
using gradients. 

Consequently, the method of the invention basically 
makes it possible : 
5 - for the uncertainty inherent in the geological 

model to be interpreted by quantifying the uncertainty 
relating to production, 

- to test whether the predicted hypotheses axe 
compatible with the data known about the reservoir from 

10 the outset. 

Finally, combining Bayesian formalism with a 
gradients method means that this method can be operated 
on the basis of more effective algorithms than those 
previously used. 

15 Other characteristics and advantages of the method of 

the invention will become clear from the following 
description of embodiments and validation tests, given by 
way of non- limitative examples, and from the attached 
drawings, in which: 

20 - Fig. 1 is a diagram of a geological model made up 

of two layers, CI, C2 , in vertical section; 

- Fig. 2 shows the a priori probability density p(kl, 
k2) ("a priori pdf") characterising the initial 
geological model as a function of permeabilities Kl, K2 

25 of the two layers CI, C2 of the model of Fig. 1; 

- Fig. 3 shows how the resemblance or likelihood 
function p(P=Pm/Kl, K2) varies, characterising the 
probability of a simulation as a function of these same 
parameters; 
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- Fig. 4 shows the variation in the a priori 
probability density p(kl, k2/P-Pm) ("a posteriori pdf" ) 
obtained by Bayesian inversion. 

- Fig. 5 shows the automatic adjustment process, 
5 based on combining Bayesian formalism with the gradients 

method, which allows the optimum of the a posteriori 
probability density to be obtained. The process is 
illustrated starting from four different starting points 
SP1, SP2 , SP3 , SP4; 
10 - Fig. 6 shows, in the parameter space, the 

geological models obtained for an optimum prediction, a 
maximum prediction and a minimum prediction; and 

- Fig. 7 shows the curve FPa representing the 
predicted variations over time of the bottom hole 

15 pressure in the production well, surrounded by the curves 

FPM and FPm corresponding to the maximum and minimum 
forecasts respectively . 

DESCRIPTION OF THE METHOD 

20 

The method of adjusting production scenarios as 
proposed by the invention consists of the following 
steps : 

1) Firstly, an initial geological model of the 
25 reservoir zone under study is constructed, incorporating 

all the information available beforehand. An inversion 
method of a known type and in particular the Bayesian 
inversion process described above can be used to carry 
out a first production history match. The adjusted model 
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can be used to obtain a mean prediction. 

2) Once this initial model has been established, one 
or several possible scenarios of developments during the 
course of production are defined. These scenarios 
correspond to hypotheses pertaining to how production 
will develop over time. Fictitious data corresponding to 
given future instants are added for each of these 
scenarios. These fictitious data can relate to a possible 
water inflow at a given time, for example, or hypotheses 
relating to oil recovery. 

3) The parameters of the initial model are adjusted 
to obtain a match with both the production history 
(observed data) and the additional data relating to tlie 
production hypotheses. A conventional inversion process 
is used for each scenario in order to obtain a new 
geological model, defined by a new set of parameters x 
and corresponding to the previous hypotheses. 

4) The probability of each scenario created is 
verified directly using an inversion formalism of tlie 
Bayesian type, for example, with the value of the a 
posteriori probability density function. 

5) The models obtained are used to produce direct 
simulations for each of these scenarios, so that new 
production forecasts can be produced. 

Several situations might arise. By considering tlie 
choice of parameters and the constraints on the 
geological model, a set of parameter values can be found 
which will enable the measured and added data to be 
reproduced. This information can be very useful to the 
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reservoir engineer. If, for example, * the aim of the 
scenario created is to test whether a premature inflow of 
water could occur in a given well and if the process of 
evaluating the scenario produces a set of physically 
5 admissible parameter values, the reservoir engineer will 

be in a position to take decisions about how the 
reservoir will be developed so that this risk can be 
prevented or he may decide to carry out new measurements 
to check whether the risk is real . 

10 Another conceivable situation is where the adjustment 

procedure can be applied to the scenario but only at the 
expense of considerable distortion to the initial 
geological model. The reservoir engineer can then be 
confident in his forecasts and be assured that the 

15 occurrence of such an event is highly unlikely. In the 

same hypothesis of preventing a premature inflow of 
water, he can regard the risk as being very low before a 
given period of time has elapsed. 

Another area where the method can be applied is for 

20 directly quantifying production uncertainties by finding 

extreme scenarios. This application will be described 
below. 

Comparison of a pessimistic scenario and an 
25 optimistic scenario as a means of assessing uncertainty 

The objective function F mentioned above in 
connection with conventional inversion techniques uses 
the information obtained from the time the reservoir was 
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placed under production to the present-. Several sets of 
parameter values may provide solutions which vrLll 
minimise the objective function. Dealing with 
uncertainties in the forecasting process from the present 
5 instant (t=tt) of the production period to a future 

instant (t=tf) may be a question of selecting from 
several possible solutions the set of parameter values 
which will give the "best" scenario and that which will 
give the "worst". The scenarios can be defined using a 

10 criterion that- is dependent on production strategy during 

the forecast period under consideration, for example 
between tt and tf . The criterion may be the total 
quantity of oil produced, for example, the total maximum 
flow rate, the instant of water inflow, the G.O.R., etc. . 

15 The best way of resolving this problem is to apply 

the method of the invention to two optimisation problems 
where the first corresponds to the "best" scenario and 
the second the "worst". 

In the first instance, a set of parameters is found 

2 0 which gives both a good match with the production history 

of the reservoir (i.e. the set of parameters which 
minimises Fm + Fx) and a maximum value for the criterion: 
Min(Fm + Fx) + Max (Criterion) . 
In the second instance, the set -of parameters is 

25 found which combines a good history match with a minimum 

value of the criterion: 

Min(Fm + Fx) + Min (Criterion) . 
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Formalism of the method: 

A future instant tf is selected and new data axe 
added to the existing production data. Each of the new 
values corresponds to a predicted production scenario. 
The predicted min/max production scenarios at this futuxe 
instant tf are determined. 

To this end, the parameter values are adjusted in 
order to obtain an optimum match of the production 
history whilst distorting the production forecast curve 
so that the minimum or maximum is attained at the future 
time t=tf. An iso-probability constraint is imposed in 
order to check the feasibility of each scenario. 

The problem is reduced to a problem of optimisation 
with non-linear constraints. The objective function used 
is 

F = F m + F x + F £ 
F ffi + F x represents the previous objective function 
incorporating the a priori data. F f can be expressed as: 
F £ = 1/2 (d f - g(x,t f )) 2 ; where 
d f represents the additional data at the time t f ; and 
g(x,t ; ) represents the production forecasts at the 
time t £ . 

F f is an additional constraint used to obtain the 
min/max scenarios. The added data d f are adjusted to 
satisfy the following two constraints: 

1) F 0 + F x < F^c, where F^ is a given criterion 
measuring the quality or the admissible standard of the 
D>*oduction history match- 



2) g(X/t f ) represents a maximum or minimum forecast. 

The first constraint is equivalent to a constraint of 
relative probability pertaining to an a posteriori 
probability density function by reference to the maximum 
probability. A ratio of a posteriori probability (ppr) is 
given : 

ppr = p (x/d-dj /p (x°°/d=dj , 
where xoo is the optimum position (parameter values 
corresponding to the maximum probability) . The admissible 
standard criterion of the history match is then: 
F^c = F m (xoo) + F x (xco) - log (ppr) 
A quadratic convergence optimisation algorithm, such 
as the Gauss -Newton algorithm known to specialists, can 
be used to minimise the objective function (F=F W + Fx + 
Ff ) , 

This algorithm can be used to carry out an automatic 
process of adjustment to the history, finding the maximum 
or the minimum of the production forecasts 

simultaneously. 

It can be of advantage to use a multi-purpose 
reservoir simulator to carry out this procedure, bringing 
the gradients method into play. By using the gradients 
method, the sensitivity of the model to the inversion 
parameters is in effect obtained directly. This method is 
well know to specialists and has been described, for 
example, by: 

Anterion F. et al . , "Use of parameter gradients for 
reservoir history matching", SPE 18433, Houston, TX, 
February 6-8, 198 9. 
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When using gradient values, the optimisation 
algorithm is very effective in dealing with the problem 
in the context of a global iterative method, 
incorporating the search for extreme scenarios and 
5 optimum matching of the production history. The method of 

the invention allows an inversion algorithm to be used 
that combines the advantages of Bayesian formalism and 
the efficiency of optimisation algorithms based on the 
use of gradients* The method enabling the Bayesian 
10 formalism to be combined with the gradients method is 

described below. 

Modified Bayesian type inversion formalism 

15 When carrying out the inversion procedures, it can be 

of advantage to use a specific inversion method derived 
from the Bayesian type of formalism, which allows a 
direct and much faster search than conventional formalism 
to find the optimum position of the a posteriori 

20 probability density function. 

The objective is achieved by defining a function to 
be optimised and by combining the Bayesian formalism with 
the gradients method. The initial information associated 
with the geological model can be incorporated in this 

25 function to be minimised. In addition, by making certain 

assumptions, it is possible to obtain an approximation of 
the a posteriori probability density around this optimal 
solution. To include this initial knowledge in the 
function to be minimised (designated as - * — ~™ 
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which is dependent on the difference between a given a 
priori average 2^ and the value x of the following 
parameters, is introduced: 

F - (F. + F x )/2, with 
5 F B « (d^-d) 7 C^fc^-d), and . 

F x = (2£a - 2c) T C x " l (2c a - x) . 
Covariance operators are introduced to model the 
Gaussian uncertainties on the measurements (matrix C m ) and 
the a priori probability density on the parameters 
10 (matrix C^O . 

The proposed formula is an extension of the 
conventional objective function of least squares. There 
is a direct link between the objective function and the 
Bayes rule. The a posteriori probability density function 
15 can be expressed in the form: 

p(x/d-dJ=Cte.exp (-F) =Cte {exp (-FJ .exp(-F x ) } , with 
P(d=d ffl /x = Cte x {exp (-F.) } 
p(x) = Cte 2 {exp (-F x )}. 
The parameters obtained for the maximum of the "a 
20 posteriori probability density function, denoted as x°°/ 

can therefore be considered to be those which minimise 
the objective function F. This correspondence is well 
known to those specialised in inversion. The modified 
procedure differs in that: 
25 - it uses the gradients method to minimise the 

objective function, which, in accordance with the 
invention, makes it possible to obtain extreme production 
scenarios with an iso-probability criterion; and in that 

- it uses the gradients method combined with Bayesian 



formalism to obtain, by linearisation, an approximation 
of the a posteriori probability density in the vicinity 
of the optimum. In practice, this method is able to meet 
satisfactorily the constraints of a posteriori iso- 
probability in the inversion algorithms. 

Once convergence has been reached in the optimisation 
procedure and the optimal values x°° have been found using 
the gradients method, assuming that the a posteriori 
probability density is locally Gaussian, it is modelled 
by an a posteriori covariance operator. Calculation of 
this covariance operator gives an approximation of the "a 
posteriori probability density" function. 

Experimental example 

The validity of the method of the invention was 
tested on a synthetic experimental case, set out below, 
which demonstrates clearly the improvement attained by 
combining inversion formalism of the Bayesian type with a 
gradient method and the possibility of quantifying 
uncertainties about the predictions by researching 
extreme scenarios . 

Having defined a synthetic geological representation, 
a numerical simulation was performed to obtain synthetic 
production data using an industrial reservoir simulator 
of the usual type . 

An initial geological representation was chosen, 
representing the information available about the 
r_e s e rvo_ij^_b_e_f _o.r_e hand_._T_he__gj eolocHLc a J moH pJ in t h i s c ase. 



is a vertical section having two layers (Fig. 1) CI, C2 , 
7.5 tn in height and 100 m thick, with constant horizontal 
permeabilities Kl, K2 of 200 mD and 100 mD respectively. 
The simulation grid is made up of ten cells for each 
layer in the horizontal direction and the total length is 
500 m. 

Water is injected, on the left-hand side of the 
reservoir, into the lower section and a total liquid flow 
rate (water + oil) is produced on the right-hand side of 
the upper layer. The same volume flow rate (200 m 3 /j or 
1250Bl/j) is imposed on the injector and on the producing 
well in order to maintain the average pressure and 
prevent the oil pressure from falling below the bubble 
point . 

The fluids that were initially present in the 
reservoir form a petroliferous reservoir containing 12% 
of irreducible water. A numerical simulation is conducted 
during a production period of 600 days. During an initial 
period, (from t = 0 to t = 300 days) , only the oil phase 
is produced. The water does not reach the producing well 
until after a period of about 300 days. The volume of 
water produced increased during the final period: 10% of 
water at t = 380 days and 77.7% of water at the end of 
the simulation (t = 600 days) . 

Synthetic measurements were produced using the 
results of this reference simulation. Only the first 
production period from t « 0 to t = 300 days) was used: 
in order to conserve a greater degree of uncertainty with 
respect to the initial characterisation of the reservoir , 



22 

it was assumed that the time at which a water inflow 
occurred was unknown. The measurement selected was down 
hole pressure in the producing well with five values 
observed between 25 and 300 days. 

An overall uncertainty range was selected to describe 
both the measurement uncertainties and the numerical 
simulation errors. These uncertainties were modelled 
using the Gaussian probability densities centred on the 
results of the reference simulation and with a=13 8 kPa 
(20psi) for the standard deviations. 

These standard deviations were implemented in the 
objective function in the diagonal terms of the 
covariance operator C m relating to the measurements. If 
one considers the first term of the objective function F m 
which quantifies the variance between the observations d^ 
(in this case the down hole pressure) and the response of 
the simulation g(x,t), the diagonal terms of C ro _1 are 
equal to l/a m 2 . 

Defining the a priori geological model 

It is assumed that the geological reference model is 
unknown and that some geological information is available 
beforehand. On a hypothetical basis, the geometry is 
considered to be well characterised, that the horizontal 
permeabilities assigned to each layer are unknown and 
that a priori information about these permeabilities is 
available but with a given margin of uncertainty. 
TKe tLw.o__h o_-r_i .2 ont a J D_e.rme.ab_i l.i.t.i.e s havie_b.e_en__s_ e 1 e c t ed. 
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as being the parameters to be subjected to inversion. The 
a priori information available is modelled by means of 
Gaussian uncertainties (Fig. 2) . The variance between the 
a priori permeabilities x^ and the parameters x is 
quantified in the second term of the objective function 
F x . The standard deviations of the parameters a x axe 
incorporated in the covariance operator C x in the same way 
as for modelling the measurement errors: the diagonal 
terms of C/ 1 are equal to l/a x 2 . The extra-diagonal terms 
are zero for this example, which means that the a priori 
data pertaining to a possible correlation between the two 
parameters were not incorporated. 

The numerical values selected a priori are set out in 
table 1 below: 
15 Table 1 : A priori geological model 

A priori value Standard deviation Value of the 

reference case 
Kl 160 mD 40 mD 200 mD 

(upper layer). 

20 K2 160 mD 40 mD 100 mD 

(lower layer) 

Comparison of different inversion' methods for thte 
test case 
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A) Automatic history matching method 

Firstly, a conventional history matching method 
combined with an automatic optimisation algorithm is 
applied. In this initial aoproach, the ceoloa^^ 1 = 
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available beforehand are used only to select the starting 
point of the process. The objective function to be 
minimised incorporates only the first term F ffi and does not 
take account of the a priori model. 
5 The initial parameters are fixed at 160 mD (a priori 

values) . The automatic matching procedure based on a 
gradient method and the conventional Gauss -Newton 
optimisation algorithm is used. In this case, the 
sensitivity of the down hole pressure in the producing 

10 well to the permeability values is obtained by the 

gradients method. 

Because synthetic measurements produced by the 
numerical simulator are used, the adapted model obtained 
in this way is very close to the reference case. Only 

15 five simulations are needed to obtain the position of the 

function minimum. The optimisation process is summarised, 
from the starting point to an optimum, in table 2 below: 

Table 2 : Automatic optimisation process for history 
matching 

2 0 Iteration 



Number 


Kl 






K2 




Fm 


0 


160 


mD 


160 mD 


0 . 


.799 


1 


172 . 


8 


mD 


133 . 6 mD 


0 . 


.407 


2 


193 . 


4 


mD 


104.9 mD 


0 . 


.042 


3 


200. 


8 


mD 


99.1 mD 


0 


. 0008 


4 


200. 


5 


mD 


99.9 mD 


0 


. 0003 
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Use of the Bayesian type inversion formalism 

By using the Bayesian inversion method, the a 
posteriori probability density is obtained, mapped over 
5 the entire parameter space. In this instance, tihe 

information available beforehand (Fig. 2) is taken into 
account in this a posteriori probability density by 
combining it with the likelihood function (Fig. 3) 
calculated for each group of parameters . 

10 The parameter space has been mapped with a large 

range of permeabilities, from 20 mD to 300 mD, with 35 
values for each parameter (the total number of numerical 
simulations will be 35 2 = 1225 operations) . It is of 
advantage of have a full map of the a posteriori 

15 likelihood function in that it enables numerous data 

about the simulation behaviour to be obtained. If we look 
at the graphical representation of the function (Fig. 4), 
the following points can be noted: 

- only one optimum is clearly identified, located at 
20 143 mD (permeability of the upper layer) and at 173 mD 

(permeability of the lower layer) . The maximum value of 
the a posteriori probability density is 1.48.10" 7 ; 

- it was possible to obtain a better characterisation 
of the geological model with a lesser degree of 

25 uncertainty with regard to the parameters (compared with 

the a priori probability density function of Fig. 2) ; 

- a correlation is established between the two 
parameters; the shape of the surface is not symmetrical 
around the oor^mrm hrr * «=; snorted on an axis 



correlation. 

On the other hand, integration of the surface makes 
it possible to calculate the a posteriori statistical 
parameters with regard to the parameters. 

Tables 4 and 5 below show the standard a posteriori 
deviations compared with the initial data: 

Table 4 ; Initial and a posteriori parameters 

Initial A posteriori 

Kl 160 mD 143 mD 

K2 160 mD 173 mD 

Table 5 : Initial and a posteriori standard deviations 

Initial A posteriori 

Kl 4 0 mD 32.8 mD 

K2 40 mD 19.0 mD 

Using Bayesian type formalism combined with the 
gradients method 

This modified inversion procedure combining, as we 
have seen, inversion formalism of the Bayesian type with 
a gradients method was applied to the experimental case. 

The optimisation algorithm is applied to find the 
optimum of the a posteriori probability density function 
directly. Four different starting points were used to 
cover the initial range of uncertainties. For each 
optimisation process, the optimum is reached with between 
4 and 7 iterations (Fig. 5) . The proposed statistical 
analysis is applied in the area around the optimum. The 
results " obtained 
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Bayesian inversion method in table 6 below: 

Table 6 : Comparative results obtained 



with 



conventional Bayesian inversion and with the modified 
inversion 



Number of operations 
Optimum (en mD) 



Bayesian 
inversion 
1225 
Kl=143 
K2=173 



Standard deviations (in mD) Kl=32.8 

K2=19 . 0 



Modified 
inversion 
20 

Kl = 143 .2 
K2=172.9 

Kl=31.4 
K2=18 . 8 



Max. of the a posteriori 
15 probability density pdf 1.48.10" 7 1.57.10" 7 

This table clearly illustrates the rate of execution 
possible when using the modified inversion method of the 
invention . 
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Using the scenario method 
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The results obtained with the experimental case are 
set out in Fig. 6. After using the automatic procedure 
described above to find the optimum of the a posteriori 
probability density and thus to obtain the match with the 
actual data (from t=0 to t=300 days) , a direct numerical 
simulation taking account of the a priori model was used 
up to t ; =600 days. This simulation gives an average 
forecast of the down hole pressure. The a priori rate of 
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probability ppr is then selected, in this case ppr=0 . 1 . 
The method of the invention is used to identify the 
extreme production scenarios, the geological model 
remaining with this range of probability. 

The constraints of iso-probability are illustrated in 
Fig. 7 (using Bayesian inversion mapping) in order to 
superpose the values corresponding to the three scenarios 
(optimum, maximum and minimum) in the parameter space. 
The location of parameters corresponding to the mini/maxi 
scenarios is close to the border defined by the rate of 
probability ppr=0.1. The numerical results are set out in 
table 7 : 

Table 7 : Extreme scenarios and comparison with the 
optimum 





Maxi . 


Optimum 


Mini . 


Pressure (psi) 


4179 


4098 


4020 


Kl (mD> 


90.2 


143 .2 


200 .8 


K2 (mD) 


217.4 


173 .0 


136.5 


The uncertainty 


concerning 


the forecast 


of the 



hole pressure at t f =600 days is given by the difference 
between the maximum forecast and the minimum forecast. 
This demonstrates that it is possible to transpose the 
uncertainties pertaining to the geological model into 
uncertainties pertaining to production forecasts. 
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CLAIMS 

1. A method for predicting developments during 
production in an underground zone containing fluids , such 
as a hydrocarbon reservoir, wherein it comprises the 

5 following steps : 

a) an initial geological model defined by a set of 
parameters x is constructed by integrating the data 
available and an inversion method is applied in order to 
determine the values of parameters that will allow the 

10 response of the model to be adjusted to a known 

production history for the said zone; 

b) one or several possible scenarios of developments 
during production are defined and, for each of these 
scenarios, observed data are added to fictitious data at 

15 given future instants corresponding to hypotheses; 

c) for each of the said scenarios, the parameters of 
the said model are adjusted to reconstruct both actual 
measurements from the production history and data added 
at future times and an inversion procedure is applied to 

20 obtain for each scenario a new set of parameters x 

characterising a modified geological model and 
corresponding to the hypotheses of the scenario; and 

d) for each of these modified geological models, 
simulations are performed to establish new production 

25 forecasts. 

2. A method as claimed in claim 1, wherein 
scenarios corresponding to extreme forecasts are 
determined so as to quantify the uncertainties of 



3- A method as claimed in claim 1 or 2, wherein 

the matching of the scenarios is adjusted in a process of 
optimisation by a combination of -a Bayesian type 
formalism and a gradients method. 

4 - A method as claimed in claim 3, wherein .a 
gradients method, is used to obtain an approximation of 
the a posteriori probability density. 

5 • A method substantially as hereinbef orre 
described with reference to the drawings . 
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